*IAT SCORING ALGORITHM FROM GREENWALD, NOSEK, & BANAJI, 2002
LAST UPDATE: 03/03/04

This macro will transform a datafile with raw trial latencies (stored as one line per response) 
for a standard format IAT (7 blocks) into a one line summary per subject of the IAT effect using
GNB's new scoring algorithm.

To use this algorithm for your SAS program, perform these two steps:  
(1) Run this script.  The macro should then be loaded into SAS's active memory.  
(2) In your SAS program, enter the following statement 

%iatCalc(libIn, libOut, indata, outdata, BLOCNAME, SESHID, TRLLATNC, TRLERROR, VERROR, VEXTREME, VSD)

and replace the dummy values with the ones that correspond to your own datafiles and variable names.  


Descriptions of what the macro expects for input, and what it will output are below.  Also, at the 
bottom of this commenting section, there is a sample SAS script that would invoke this macro (assuming 
that you have already run the macro to get it into SAS's active memory).

   The macro expects the following types of libraries, datafiles, and data:
      libIn     = input SAS library name
      libOut    = output SAS library name
      indata    = filename of input SAS dataset in the input SAS library
      outdata   = filename for the output SAS dataset in the output SAS library
      BLOCNAME  = variable name for block identifier in the indata file: alphanumeric indication of the 
                  four trial blocks ('B3', 'B4', 'B6', 'B7' are critical blocks corresponding 
                  to B3, B4, B6, and B7 from GNB, 2002).  At present the macro requires that the variable 
                  passed here uses the names 'B3', 'B4', 'B6', 'B7' to refer to the B3-B7 blocks
	  SESHID    = variable name for unique subject identifier in the indata file
	  TRLLATNC  = variable name for latency of response for trial in the indata file
	  TRLERROR  = variable name for error coding: 0 if initial response was correct, 1 if initial response 
                  was incorrect in the indata file
      VERROR    = value: if '1' the algorithm will use error trial latencies, if '2' the algorithm will replace 
                  error trial latencies with blockmean+600
      VEXTREME  = value: if '1' the algorithm provides no treatment of extreme values, if '2' the algorithm
                  will delete trials <400ms
      VSTD      = value: if '1' the block standard deviation is performed including corrected error trials,
                  if '2' the block standard deviation is performed on correct responses only (1 is standard)

   The macro will output the following variables to the outdata file:
      SESHID  = unique subject identifier
      SUBEXCL = 0 for inclusion data, 1 for excluded data, 2 for incomplete data
      MB3     = mean of trial latencies for B3
      MB4     = mean of trial latencies for B4
      MB6     = mean of trial latencies for B6
      MB7     = mean of trial latencies for B7
      CS1     = standard deviation for B3 and B6 trials combined (correct trials only)
      CS2     = standard deviation for B4 and B7 trials combined (correct trials only)
      AS1     = standard deviation for B3 and B6 trials combined (all trials)
      AS2     = standard deviation for B4 and B7 trials combined (all trials)
      EB3     = percent errors of trials for B3
      EB4     = percent errors of trials for B4
      EB6     = percent errors of trials for B6
      EB7     = percent errors of trials for B7
      FB3     = percent fast responses of trials for B3
      FB4     = percent fast responses of trials for B4
      FB6     = percent fast responses of trials for B6
      FB7     = percent fast responses of trials for B7
      DIFF1   = MB6 - MB3
      DIFF2   = MB7 - MB4
      IAT1    = DIFF1/STD1
      IAT2    = DIFF2/STD2
	  IAT     = mean of IAT1 and IAT2
;

/* SAMPLE SAS SCRIPT THAT USES THIS MACRO ASSUMING THAT IT IS IN ACTIVE SAS MEMORY
This sample script opens two libraries called 'web' and 'outdata'. In the 'web' library a SAS file 
named 'iatrace350000' already exists and will be pulled into active memory by the macro.  When the macro 
finishes it will output the result data into a file in the library 'outdata' called 'Ciatrace350000'. The 
information coding the name of the block is in the variable BLOCK (which is passed to 'BLOCNAME' for
the algorithm). The information coding subject ID is SUB (passed to 'SESHID').  The information coding
the latency of each trial is in LATENCY (passed to TRLLATNC).  The information coding whether an error
was made on each trial is in ERROR (passed to TRLERROR).  And, the setting '2' for VERROR means that 
error latencies will be replaced with the blockmean+600.  The setting '2' for VEXTREME means that responses
faster than 400ms will be deleted before calculating the IAT effect. The setting '1' for VSTD means that 
the standard deviations will be calculated after recombining error trials so the range of D values is bounded 
between (-2,2).  When the algorithm finishes and saves the output, the proc means statement will calculate a 
mean of the IAT effects across subjects (after excluding 'bad' subjects).

libname web 'H:\raceatt\';
libname outdata 'H:\raceattclean\';
%iatCalc(web, outdata, iatrace350000, Ciatrace350000, BLOCK, SUB, LATENCY, ERROR, 2, 2,1);
proc means (data=outdata.Ciatrace350000);var IAT;where SUBEXCL=0;run;
*/

/*MACROS for NEW IMPLICIT ASSOCIATION TEST SCORING ALGORITHM (Greenwald, Nosek, & Banaji, 2002)*/
%macro iatCalc(libIn, libOut, indata, outdata, BLOCNAME, SESHID, TRLLATNC, TRLERROR, VERROR, VEXTREME, VSTD);
%let divide= /; %let multiply= *; %let add = +; %let subtract = -;
%iatAlgorithm(libIn=&libIn, libOut=&libOut, indata=&indata, outdata=&outdata, BLOCNAME=&BLOCNAME, 
              SESHID=&SESHID, TRLLATNC=&TRLLATNC, TRLERROR=&TRLERROR, VERROR=&VERROR, VEXTREME=&VEXTREME, VSTD=&VSTD);run;
%mend iatCalc;

%macro iatAlgorithm(libIn, libOut, indata, outdata, BLOCNAME, SESHID, TRLLATNC, TRLERROR, VERROR, VEXTREME, VSTD);
data IAT; set &libIn..&indata; 
%*PRELIMINARY STEPS FOR HANDLING WEBDATA FORMATS;
    keep &BLOCNAME &SESHID &TRLLATNC &TRLERROR;
proc sort data=iat; by &SESHID &BLOCNAME;	
	%*options nonotes; %*suppress all Notes to log;	

data IAT; set IAT;

%*STEP 1: Include data from B3, B4, B6, B7;
   if &BLOCNAME in ('B3', 'B4', 'B6', 'B7') then ;
   else delete;

%*STEP 2a: Eliminate trial latencies > 10,000ms;
   if &TRLLATNC > 10000 then delete;

%*STEP 2b: Eliminate subjects for whom more than 10% of trials have latencies < 300ms;
   else if &TRLLATNC < 0 then delete; %*for miscoded data in datafile indicating negative response times;
   else if -1 < &TRLLATNC < 300 then FAST = 1;
   else FAST = 0;

data FASTDATA; set IAT; keep &SESHID &BLOCNAME FAST;
proc sort; by &SESHID &BLOCNAME;
proc means data=IAT noprint; by &SESHID &BLOCNAME; var FAST; output out=means mean=MEAN;
proc transpose data=means prefix=F name=name out=FASTMEAN; by &SESHID; id &BLOCNAME;
data FASTMEAN; set FASTMEAN; where name='MEAN'; FASTM = mean(FB3, FB4, FB6, FB7);
   if FASTM > .10 then SUBEXCL = 1; else SUBEXCL = 0; 
   %*SUBEXCL = 0 (include data), 1 (exclude data - too many fast responses), 2 (exclude data - missing data);
   %*The SUBEXCL variable needed to be reintroduced to the final dataset in STEP 12;

%*STEP 3: Use all trials;
   %*in the conventional algorithm the first two trials of each block would be dropped here;

%*STEP 4: No extreme value treatment <or> delete trial with latencies <400ms;
data IAT; set IAT; 
    %*if &VEXTREME = 1 then do nothing here;
    if &VEXTREME = 2 then do; if &TRLLATNC < 400 then delete; end;

proc sort data=IAT; by &SESHID &BLOCNAME;

%*STEP 5: Compute mean of correct latencies for each block;
data CORR; set IAT;
   %*if &VERROR is 1 then means and SDs will be calculated for the entire set of latencies;
   if &VERROR = 2 then do; if &TRLERROR NE 0 then delete; end;
   keep &SESHID &BLOCNAME &TRLLATNC &TRLERROR;
proc means data=CORR noprint; by &SESHID &BLOCNAME; var &TRLLATNC; output out=means mean=MEAN;
proc transpose data=means prefix=M name=name out=CORRMEAN; by &SESHID; id &BLOCNAME;
data CORRMEAN; set CORRMEAN; where name='MEAN';

%*STEP 6a: Compute pooled SD for B3 & B6, and separately for B4 & B7 for correct trials only;
data SD; set CORR;
   if &TRLERROR NE 0 then delete;  *drop error trials;
   if &BLOCNAME in ('B3', 'B6') then TD = '1';
   else if &BLOCNAME in ('B4', 'B7') then TD = '2';
   else delete;
   drop &BLOCNAME;
proc sort data=SD; by &SESHID TD;
proc means data=SD noprint; by &SESHID TD; var &TRLLATNC; output out=means std=STD;
proc transpose data=means prefix=CS name=name out=CORRSTD2; by &SESHID; id TD;
data CORRSTD2; set CORRSTD2; where name='STD';

%*STEP 6b: Compute pooled SD for B3 & B6, and separately for B4 & B7 including error trials;
data SD; merge IAT CORRMEAN; by &SESHID;
   if &TRLERROR < 0 then delete;
   else if &TRLERROR > 1 then delete; %*get rid of coding errors;
   else if &TRLERROR = 1 and &VERROR = 2 then do;
      if &BLOCNAME in ('B3') then &TRLLATNC = MB3 + 600;
	  else if &BLOCNAME in ('B4') then &TRLLATNC = MB4 + 600;
	  else if &BLOCNAME in ('B6') then &TRLLATNC = MB6 + 600;
	  else if &BLOCNAME in ('B7') then &TRLLATNC = MB7 + 600;
   end;
   if &BLOCNAME in ('B3', 'B6') then TD = '1';
   else if &BLOCNAME in ('B4', 'B7') then TD = '2';
   else delete;
   drop &BLOCNAME;
proc sort data=SD; by &SESHID TD;
proc means data=SD noprint; by &SESHID TD; var &TRLLATNC; output out=means std=STD;
proc transpose data=means prefix=AS name=name out=CORRSTD1; by &SESHID; id TD;
data CORRSTD1; set CORRSTD1; where name='STD';

%*STEP 7: Replace error latencies with block mean + 600ms 
<or> use latency from stimulus onset to correct response (when correct response is required);
data ERR; set IAT;
   keep &SESHID &BLOCNAME &TRLERROR;
   if &TRLERROR < 0 then delete;
   else if &TRLERROR > 1 then delete; %*get rid of coding errors;
proc means data=ERR noprint; by &SESHID &BLOCNAME; var &TRLERROR; output out=means mean=MEAN;
proc transpose data=means prefix=E name=name out=ERRMEAN; by &SESHID; id &BLOCNAME;
data ERRMEAN; set ERRMEAN; where name='MEAN';

%*STEP 7 continued: combining data;
data COMBINE; merge CORRMEAN CORRSTD1 CORRSTD2 ERRMEAN FASTMEAN; by &SESHID; %*combining datasets for calculating final means;
   if &VERROR=2 then do;
      array BLOCKMeans(*) MB: ;
      array BLOCKCstds(*) CS: ;
	  array BLOCKAstds(*) AS: ;
      array BLOCKErrs(*) EB: ;
      do i=1 to dim(BLOCKMeans); %*for each of the four blocks replace error trials with mean + 600ms;
	     BLOCKMeans{i} = (1-BLOCKErrs{i})*BLOCKMeans{i} + (BLOCKErrs{i})*(BLOCKMeans{i}+600);			
      end;
   end;

%*STEP 8: No transformation of latencies;
    %*in the conventional algorithm, raw latencies would be log transformed prior to the transposing in the current format;

%*STEP 9: Average latencies for each of the four blocks;
    %*this step was already accomplished in the do loop above;

%*STEP 10: Compute two differences B6-B3 and B7-B4 (does not account for pairing order);
   DIFF1 = MB6 - MB3;
   DIFF2 = MB7 - MB4;

%*STEP 11: Divide each difference by associated pooled SD from STEP 6a or 6b;
   If &VSTD = 2 then do;
      IAT1 = DIFF1/CS1;
      IAT2 = DIFF2/CS2;
   end;
   else do; %*IF VSTD = 1 also set as default;
      IAT1 = DIFF1/AS1;
	  IAT2 = DIFF2/AS2;
   end; 

%*STEP 12: Average quotients from STEP 11;
   IAT = mean(IAT1, IAT2);

%*if there is missing data, do mark data as excluded (SUBEXCL=2);
   do i=1 to dim(BLOCKMeans);
      if BLOCKMeans{i} = . then SUBEXCL = 2;
   end;
data &libout..&outdata (drop=i name);
	set COMBINE;
run;
%mend iatAlgorithm;
*END OF ALGORITHM;

libname study 'C:\primary\Data\Dissertation.2001\Diss.1\';
data fuzzy1;
  informat Code $18.;informat Bname $18.;infile 'C:\primary\Data\Dissertation.2001\Diss.1\diss.study1b.dat' delimiter='09'x;
  input date    time $ build sub     code $  block   trial   lat     resp    corr    stim1 bname $;
data fuzzy2;
  informat Code $18.;informat Bname $18.;infile 'C:\primary\Data\Dissertation.2001\Diss.1\diss.study1d.dat' delimiter='09'x;
  input date    time $ build sub   code $  block   trial   lat     resp    corr   stim1  bname $;
data fuzzy3;
  informat Code $18.;informat Bname $18.;infile 'C:\primary\Data\Dissertation.2001\Diss.1\diss.study1e.dat' delimiter='09'x;
  input date    time $ build sub   code $  block   trial   lat     resp    corr   stim1  bname $;
data fuzzy4;
  informat Code $18.;informat Bname $18.;infile 'C:\primary\Data\Dissertation.2001\Diss.1\diss.study1f.dat' delimiter='09'x;
  input date    time $ build sub   code $  block   trial   lat     resp    corr   stim1  bname $;
data fuzzy5;
  informat Code $18.;informat Bname $18.;infile 'C:\primary\Data\Dissertation.2001\Diss.1\diss.study1c.dat' delimiter='09'x;
  input date    time $ build sub   code $  block   trial   lat     resp    corr   stim1  bname $;
data fuzzy6;
  informat Code $18.;informat Bname $18.;infile 'C:\primary\Data\Dissertation.2001\Diss.1\diss.study1a.dat' delimiter='09'x;
  input date    time $ build sub     code $  block   trial   lat     resp    corr    stim1 bname $;
data fuzzy6;
  informat Code $18.;informat Bname $18.;infile 'C:\primary\Data\Dissertation.2001\Diss.1\diss.study1a.dat' delimiter='09'x;
  input date    time $ build sub     code $  block   trial   lat     resp    corr    stim1 bname $;
data fuzzy7;
  informat Code $18.;informat Bname $18.;infile 'C:\primary\Data\Dissertation.2001\Diss.1\diss.1.dat' delimiter='09'x;
  input date    time $ build sub     code $  block   trial   lat     resp    corr    stim1 bname $;
data comb;
  set fuzzy1 fuzzy2 fuzzy3 fuzzy4 fuzzy5 fuzzy6 fuzzy7; drop build time date; if resp < 15 then resp = resp - 1;
  if resp = 79 then resp = 1; if resp = 80 then resp = 2; if resp = 81 then resp = 3;
  if resp = 75 then resp = 4; if resp = 76 then resp = 5; if resp = 77 then resp = 6;
  if resp = 71 then resp = 7; if resp = 72 then resp = 8; if resp = 73 then resp = 9;
  if resp = 16 then delete; /*control quit*/
  if sub < 5 then delete;
  err = 1 - corr;

  if bname in ('bcode', 'boggle1', 'boggle2', 'boggle3', 'interview', 'badgood', 'flowersinsects', 'insectsflowers') then delete;
  if bname in ('badgood', 'blackbad', 'blackgood', 'democratbad', 'democratgood',
               'flowersbad', 'flowersgood', 'flowersinsects', 'gaybad', 'gaygood',
               'insectsflowers', 'sciencebad', 'sciencegood') then task = 1; else task = 0;

/*proc means;class bname; var lat;run;
proc means; class code;var lat;run;
proc means; class sub; var corr;run;
proc univariate;var lat;run;
proc means mean std; class sub bname; var corr;run;
proc means; class block;var lat;run;*/

proc sort; by sub;
data study.raw; set comb;

data study.iFI; set study.raw;
  if bname = 'flowersgood' and trial < 22 then blockn = 'B3';
  else if bname = 'flowersgood' and trial > 22 then blockn = 'B4';
  else if bname = 'flowersbad' and trial < 22 then blockn = 'B6';
  else if bname = 'flowersbad' and trial > 22 then blockn = 'B7';
  else delete;
%iatCalc(study, study, iFI, iFIclean, blockn, SUB, LAT, err, 2, 2, 1);
data study.iFIclean; set study.iFIclean;
  FIiat = IAT; FIiat1 = IAT1; FIiat2 = IAT2;
  FIiatexcl = subexcl; FItradIAT = DIFF2;
  FImb3 = mb3; FImb4 = mb4; FImb6 = mb6; FImb7 = mb7; 
  FIeb3 = eb3; FIeb4 = eb4; FIeb6 = eb6; FIeb7 = eb7; 
  if FIeb3>.39 | FIeb4>.39 | FIeb6>.39 | FIeb7>.39 then FIiatexcl=1;
  if FIiatexcl = 1 then do; FIiat = .; FIiat1 = .; FIiat2 = .; FItradIAT = .; end;
  drop IAT IAT1 IAT2 subexcl mb3 mb4 mb6 mb7 CS1 CS2 AS1 AS2 EB3 EB4 EB6 EB7 FB3 FB4 FB6 FB7 FASTM DIFF1 DIFF2;
*proc means;run;

data study.iDR; set study.raw;
  if bname = 'democratgood' and trial < 22 then blockn = 'B3';
  else if bname = 'democratgood' and trial > 22 then blockn = 'B4';
  else if bname = 'democratbad' and trial < 22 then blockn = 'B6';
  else if bname = 'democratbad' and trial > 22 then blockn = 'B7';
  else delete;
%iatCalc(study, study, iDR, iDRclean, blockn, SUB, LAT, err, 2, 2, 1);
data study.iDRclean; set study.iDRclean;
  DRiat = IAT; DRiat1 = IAT1; DRiat2 = IAT2;
  DRiatexcl = subexcl; DRtradIAT = DIFF2;
  DRmb3 = mb3; DRmb4 = mb4; DRmb6 = mb6; DRmb7 = mb7; 
  DReb3 = eb3; DReb4 = eb4; DReb6 = eb6; DReb7 = eb7; 
  if DReb3>.39 | DReb4>.39 | DReb6>.39 | DReb7>.39 then DRiatexcl=1;
  if DRiatexcl = 1 then do; DRiat = .; DRiat1 = .; DRiat2 = .; DRtradIAT = .; end;
  drop IAT IAT1 IAT2 subexcl mb3 mb4 mb6 mb7 CS1 CS2 AS1 AS2 EB3 EB4 EB6 EB7 FB3 FB4 FB6 FB7 FASTM DIFF1 DIFF2;
*proc means;run;

data study.iGS; set study.raw;
  if bname = 'gaybad' and trial < 22 then blockn = 'B3';
  else if bname = 'gaybad' and trial > 22 then blockn = 'B4';
  else if bname = 'gaygood' and trial < 22 then blockn = 'B6';
  else if bname = 'gaygood' and trial > 22 then blockn = 'B7';
  else delete;
%iatCalc(study, study, iGS, iGSclean, blockn, SUB, LAT, err, 2, 2, 1);
data study.iGSclean; set study.iGSclean;
  GSiat = IAT; GSiat1 = IAT1; GSiat2 = IAT2;
  GSiatexcl = subexcl; GStradIAT = DIFF2;
  GSmb3 = mb3; GSmb4 = mb4; GSmb6 = mb6; GSmb7 = mb7; 
  GSeb3 = eb3; GSeb4 = eb4; GSeb6 = eb6; GSeb7 = eb7; 
  if GSeb3>.39 | GSeb4>.39 | GSeb6>.39 | GSeb7>.39 then GSiatexcl=1;
  if GSiatexcl = 1 then do; GSiat = .; GSiat1 = .; GSiat2 = .; GStradIAT = .; end;
  drop IAT IAT1 IAT2 subexcl mb3 mb4 mb6 mb7 CS1 CS2 AS1 AS2 EB3 EB4 EB6 EB7 FB3 FB4 FB6 FB7 FASTM DIFF1 DIFF2;
*proc means;run;

data study.iBW; set study.raw;
  if bname = 'blackbad' and trial < 22 then blockn = 'B3';
  else if bname = 'blackbad' and trial > 22 then blockn = 'B4';
  else if bname = 'blackgood' and trial < 22 then blockn = 'B6';
  else if bname = 'blackgood' and trial > 22 then blockn = 'B7';
  else delete;
%iatCalc(study, study, iBW, iBWclean, blockn, SUB, LAT, err, 2, 2, 1);
data study.iBWclean; set study.iBWclean;
  BWiat = IAT; BWiat1 = IAT1; BWiat2 = IAT2;
  BWiatexcl = subexcl; BWtradIAT = DIFF2;
  BWmb3 = mb3; BWmb4 = mb4; BWmb6 = mb6; BWmb7 = mb7; 
  BWeb3 = eb3; BWeb4 = eb4; BWeb6 = eb6; BWeb7 = eb7; 
  if BWeb3>.39 | BWeb4>.39 | BWeb6>.39 | BWeb7>.39 then BWiatexcl=1;
  if BWiatexcl = 1 then do; BWiat = .; BWiat1 = .; BWiat2 = .; BWtradIAT = .; end;
  drop IAT IAT1 IAT2 subexcl mb3 mb4 mb6 mb7 CS1 CS2 AS1 AS2 EB3 EB4 EB6 EB7 FB3 FB4 FB6 FB7 FASTM DIFF1 DIFF2;
*proc means;run;

data study.iSH; set study.raw;
  if bname = 'sciencebad' and trial < 22 then blockn = 'B3';
  else if bname = 'sciencebad' and trial > 22 then blockn = 'B4';
  else if bname = 'sciencegood' and trial < 22 then blockn = 'B6';
  else if bname = 'sciencegood' and trial > 22 then blockn = 'B7';
  else delete;
%iatCalc(study, study, iSH, iSHclean, blockn, SUB, LAT, err, 2, 2, 1);
data study.iSHclean; set study.iSHclean;
  SHiat = IAT; SHiat1 = IAT1; SHiat2 = IAT2;
  SHiatexcl = subexcl; SHtradIAT = DIFF2;
  SHmb3 = mb3; SHmb4 = mb4; SHmb6 = mb6; SHmb7 = mb7; 
  SHeb3 = eb3; SHeb4 = eb4; SHeb6 = eb6; SHeb7 = eb7; 
  if SHeb3>.39 | SHeb4>.39 | SHeb6>.39 | SHeb7>.39 then SHiatexcl=1;
  if SHiatexcl = 1 then do; SHiat = .; SHiat1 = .; SHiat2 = .; SHtradIAT = .; end;
  drop IAT IAT1 IAT2 subexcl mb3 mb4 mb6 mb7 CS1 CS2 AS1 AS2 EB3 EB4 EB6 EB7 FB3 FB4 FB6 FB7 FASTM DIFF1 DIFF2;
*proc means;run;


data exptemp; set study.raw; if task = 1 then delete; 
proc sort; by sub bname stim1;run;
data exptemp; set exptemp; drop stim1 code block corr task trial lat;
proc transpose data=exptemp prefix=q name=name out=compexp; by sub; run;
data compexp; set compexp; if name = 'err' then delete;

data compexp; set compexp;

if mod(sub,16) in (1,2,3,4,5,6,7,8) then do;
  b_aa1 = 10-q1; b_aa2 = q2; b_aa3 = 10-q3; b_aa4 = q4; m_aa1 = q5; m_aa2 = q6; m_aa3 = q7; m_aa4 = q8;
  b_dm1 = q9; b_dm2 = q10; b_dm3 = q11; b_dm4 = 10-q12; b_dm5 = 10-q13; b_dm6 = 10-q14; b_dm7 = q15; b_dm8 = q16;
  m_dm1 = q17; m_dm2 = q18; m_dm3 = q19; m_dm4 = q20; b_gy1 = 10-q21; b_gy2 = q22; b_gy3 = 10-q23; b_gy4 = 10-q24; b_gy5 = q25;
  m_gy1 = q26; m_gy2 = q27; m_gy3 = q28; m_gy4 = q29; b_hm1 = q30; b_hm2 = q31; b_hm3 = q32; b_hm4 = q33;
  m_hm1 = q34; m_hm2 = q35; m_hm3 = q36; m_hm4 = q37; sc_pu = q38; hm_pu = q39; gy_pu = q40; sp_pu = q41;
  aa_pu = q42; wh_pu = q43; dm_pu = q44; rp_pu = q45; sc_pn = q46; hm_pn = q47; gy_pn = q48; sp_pn = q49;
  aa_pn = q50; wh_pn = q51; dm_pn = q52; rp_pn = q53; m_rp1 = q54; m_rp2 = q55; m_rp3 = q56; m_rp4 = q57;
  b_sc1 = q58; b_sc2 = q59; b_sc3 = q60; b_sc4 = 10-q61; m_sc1 = q62; m_sc2 = q63; m_sc3 = q64; m_sc4 = q65;
  m_sp1 = q66; m_sp2 = q67; m_sp3 = q68; m_sp4 = q69; m_wh1 = q70; m_wh2 = q71; m_wh3 = q72; m_wh4 = q73;
  sc_gb = .; hm_gb = .; gy_gb = .; sp_gb = .; aa_gb = .; wh_gb = .; dm_gb = .; rp_gb = .;
  sc_wc = .; hm_wc = .; gy_wc = .; sp_wc = .; aa_wc = .; wh_wc = .; dm_wc = .; rp_wc = .;end;
else do;  /*subs 9-16*/
  b_aa1 = 10-q1; b_aa2 = 10-q2; b_aa3 = 10-q3; b_aa4 = q4; m_aa1 = q5; m_aa2 = q6; m_aa3 = q7; m_aa4 = q8;
  b_dm1 = 10-q9; b_dm2 = q10; b_dm3 = 10-q11; b_dm4 = 10-q12; b_dm5 = q13; b_dm6 = q14; b_dm7 = 10-q15; b_dm8 = q16;
  m_dm1 = q17; m_dm2 = q18; m_dm3 = q19; m_dm4 = q20; b_gy1 = q21; b_gy2 = 10-q22; b_gy3 = 10-q23; b_gy4 = 10-q24; b_gy5 = q25;
  m_gy1 = q26; m_gy2 = q27; m_gy3 = q28; m_gy4 = q29; sc_gb = q30; hm_gb = q31; gy_gb = q32; sp_gb = q33;
  aa_gb = q34; wh_gb = q35; dm_gb = q36; rp_gb = q37; b_hm1 = 10-q38; b_hm2 = q39; b_hm3 = 10-q40; b_hm4 = q41;
  m_hm1 = q42; m_hm2 = q43; m_hm3 = q44; m_hm4 = q45; m_rp1 = q46; m_rp2 = q47; m_rp3 = q48; m_rp4 = q49;
  b_sc1 = 10-q50; b_sc2 = q51; b_sc3 = q52; b_sc4 = 10-q53; m_sc1 = q54; m_sc2 = q55; m_sc3 = q56; m_sc4 = q57;
  m_sp1 = q58; m_sp2 = q59; m_sp3 = q60; m_sp4 = q61; sc_wc = q62; hm_wc = q63; gy_wc = q64; sp_wc = q65;
  aa_wc = q66; wh_wc = q67; dm_wc = q68; rp_wc = q69; m_wh1 = q70; m_wh2 = q71; m_wh3 = q72; m_wh4 = q73;
  sc_pu = .; hm_pu = .; gy_pu = .; sp_pu = .; aa_pu = .; wh_pu = .; dm_pu = .; rp_pu = .;
  sc_pn = .; hm_pn = .; gy_pn = .; sp_pn = .; aa_pn = .; wh_pn = .; dm_pn = .; rp_pn = .; end;
  drop q1-q73;

if sub=16 then hm_gb=.;if sub=51 then hm_pn=.;if sub=15 then m_dm4=.;if sub=43 then b_dm3=.;

sc_pvt = mean(sc_pu, sc_pn, sc_gb, sc_wc)-5; hm_pvt = mean(hm_pu, hm_pn, hm_gb, hm_wc)-5;
gy_pvt = mean(gy_pu, gy_pn, gy_gb, gy_wc)-5; sp_pvt = mean(sp_pu, sp_pn, sp_gb, sp_wc)-5;
aa_pvt = mean(aa_pu, aa_pn, aa_gb, aa_wc)-5; wh_pvt = mean(wh_pu, wh_pn, wh_gb, wh_wc)-5;
dm_pvt = mean(dm_pu, dm_pn, dm_gb, dm_wc)-5; rp_pvt = mean(rp_pu, rp_pn, rp_gb, rp_wc)-5;
sc_pvt1 = mean(sc_pu, sc_wc)-5; hm_pvt1 = mean(hm_pu, hm_wc)-5;
gy_pvt1 = mean(gy_pu, gy_wc)-5; sp_pvt1 = mean(sp_pu, sp_wc)-5;
aa_pvt1 = mean(aa_pu, aa_wc)-5; wh_pvt1 = mean(wh_pu, wh_wc)-5;
dm_pvt1 = mean(dm_pu, dm_wc)-5; rp_pvt1 = mean(rp_pu, rp_wc)-5;
sc_pvt2 = mean(sc_pn, sc_gb)-5; hm_pvt2 = mean(hm_pn, hm_gb)-5;
gy_pvt2 = mean(gy_pn, gy_gb)-5; sp_pvt2 = mean(sp_pn, sp_gb)-5;
aa_pvt2 = mean(aa_pn, aa_gb)-5; wh_pvt2 = mean(wh_pn, wh_gb)-5;
dm_pvt2 = mean(dm_pn, dm_gb)-5; rp_pvt2 = mean(rp_pn, rp_gb)-5;

data compexp; set compexp;
sc_i = mean(m_sc2, m_sc4); hm_i = mean(m_hm2, m_hm4);gy_i = mean(m_gy2, m_gy4); sp_i = mean(m_sp2, m_sp4);
aa_i = mean(m_aa2, m_aa4); wh_i = mean(m_wh2, m_wh4);dm_i = mean(m_dm2, m_dm4); rp_i = mean(m_rp2, m_rp4);
sc_e = mean(m_sc1, m_sc3); hm_e = mean(m_hm1, m_hm3);gy_e = mean(m_gy1, m_gy3); sp_e = mean(m_sp1, m_sp3);
aa_e = mean(m_aa1, m_aa3); wh_e = mean(m_wh1, m_wh3);dm_e = mean(m_dm1, m_dm3); rp_e = mean(m_rp1, m_rp3);
baa_pvt = mean(b_aa1, b_aa2, b_aa3, b_aa4)-5;
bdm_pvt = mean(b_dm1, b_dm2, b_dm3, b_dm4, b_dm5, b_dm6, b_dm7, b_dm8)-5;
bhm_pvt = mean(b_hm1, b_hm2, b_hm3, b_hm4)-5;
bsc_pvt = mean(b_sc1, b_sc2, b_sc3, b_sc4)-5;
bgy_pvt = mean(b_gy1, b_gy2, b_gy3, b_gy4, b_gy5)-5;

DRexp = dm_pvt - rp_pvt;SHexp = hm_pvt - sc_pvt;
GSexp = sp_pvt - gy_pvt;BWexp = wh_pvt - aa_pvt;
DRexp1 = dm_pvt1 - rp_pvt1;SHexp1 = hm_pvt1 - sc_pvt1;
GSexp1 = sp_pvt1 - gy_pvt1;BWexp1 = wh_pvt1 - aa_pvt1;
DRexp2 = dm_pvt2 - rp_pvt2;SHexp2 = hm_pvt2 - sc_pvt2;
GSexp2 = sp_pvt2 - gy_pvt2;BWexp2 = wh_pvt2 - aa_pvt2;

drop sc_pu sc_pn sc_gb sc_wc hm_pu hm_pn hm_gb hm_wc gy_pu gy_pn gy_gb gy_wc
sp_pu sp_pn sp_gb sp_wc aa_pu aa_pn aa_gb aa_wc wh_pu wh_pn wh_gb wh_wc
dm_pu dm_pn dm_gb dm_wc rp_pu rp_pn rp_gb rp_wc;

data study.total; merge compexp study.iFIclean study.iDRclean study.iGSclean study.iBWclean study.iSHclean; by sub;

*implicit multi-trait;
proc corr; var FIiat DRiat SHiat GSiat BWiat;run;

*implicit internal consistency;
proc corr; var FIiat1 DRiat1 SHiat1 GSiat1 BWiat1;
           with FIiat2 DRiat2 SHiat2 GSiat2 BWiat2;run;

*explicit multi-trait;
proc corr; var DRexp SHexp GSexp BWexp;run;

*explicit internal consistency;
proc corr; var DRexp1 SHexp1 GSexp1 BWexp1; with DRexp2 SHexp2 GSexp2 BWexp2;run;

*implicit-explicit comparison;
proc corr; var FIiat DRiat SHiat GSiat BWiat;
           with DRexp SHexp GSexp BWexp;run;



/*proc print; var sub SHiatexcl BWiatexcl GSiatexcl DRiatexcl;run;
proc means;run;*/

/*proc means; var sc_pvt hm_pvt gy_pvt sp_pvt aa_pvt wh_pvt dm_pvt rp_pvt baa_pvt bdm_pvt bhm_pvt bsc_pvt bgy_pvt;run;
proc corr; var sc_pvt hm_pvt gy_pvt sp_pvt aa_pvt wh_pvt dm_pvt rp_pvt;run;
proc corr; var baa_pvt bdm_pvt bhm_pvt bsc_pvt bgy_pvt;run;
proc corr; var sc_pvt hm_pvt gy_pvt sp_pvt aa_pvt wh_pvt dm_pvt rp_pvt;
           with baa_pvt bdm_pvt bhm_pvt bsc_pvt bgy_pvt;run;*/